{
 "cells": [
  {
   "cell_type": "markdown",
   "id": "94fb5a32",
   "metadata": {},
   "source": [
    "# 第 1 章 线性方程组的求解\n",
    "\n",
    "线性代数(Linear Algebra)顾名思义, 是研究线性运算的抽象. 这里线性运算, 表示**加法**和**数乘**, 而*代数*表示抽象和一般化. 所以这门课可以认为是在一般和抽象情形下讨论只有加法和数乘的数学模型及其求解. 这里的一般和抽象其实还不能算非常一般和抽象, 比如我们只讨论有限维(维数是 $n$)的情形. 但这种级别的抽象讨论, 已经足以使我们在面对绝大多数科学和工程的计算模型时游刃有余. 顺便也能提升我们的思维强度和深度, 或者说, 通过后天训练使你更加聪明(智力+1). \n",
    "\n",
    "尽管我们强调抽象, 但抽象不是凭空建立的. 我们要从一个具体的问题和模型入手, 然后讨论其一般形式. 也就是有一个从具体问题到一般模型和理论, 然后再到能解更多的具体问题这样的反复循环, 从而使我们不断变强. 在本章, 我们从最具体的实数域线性方程组入手, 开始我们讨论."
   ]
  },
  {
   "cell_type": "markdown",
   "id": "b5d39bd9",
   "metadata": {},
   "source": [
    "## 1.1 线性方程组的形式及其相关概念\n",
    "先把具体做到极端, 来一个我们在初中就已经熟悉二元一次方程组:\n",
    "$$\n",
    "\\left\\{\n",
    "\\begin{array}{rcl}\n",
    "3 x + 4 y &=& 2 \\\\\n",
    "x - 2 y &=& -3.\n",
    "\\end{array}\n",
    "\\right.\n",
    "$$\n",
    "然后老师曾经介绍过我们可以用消元和代入求解这类问题. 所以这个问题我们都是熟悉的. "
   ]
  },
  {
   "cell_type": "markdown",
   "id": "97ec1f37",
   "metadata": {},
   "source": [
    "对于这个二元一次方程组, 其几何意义为:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "91b5c8d2",
   "metadata": {},
   "outputs": [],
   "source": [
    "import numpy as np\n",
    "from matplotlib import pyplot as plt\n",
    "eps = np.finfo(np.float64).eps"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "63ed5f3a",
   "metadata": {
    "scrolled": true
   },
   "outputs": [],
   "source": [
    "x = np.linspace(-2,2,100)\n",
    "fig, ax = plt.subplots()\n",
    "line1, = plt.plot(x, (2 - 3*x) / 4, label = r'$3 x + 4 y = 2$')\n",
    "line2, = plt.plot(x, (-3 - x) / -2, label = r'$x - 2 y = -3$')\n",
    "plt.plot([-0.8], [1.1], 'o')\n",
    "ax.legend(handles=[line1, line2])\n",
    "ax.annotate(r'solution: $(-0.8, 1.1)$', xy=(-0.7, 1.1), xytext=(-0.0, 1.0),\n",
    "            arrowprops=dict(arrowstyle=\"->\", connectionstyle=\"arc3\"))"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "9da1843b",
   "metadata": {},
   "source": [
    "注意这个几何意义尽管直观, 对我们求解并没有什么用. 真正的求解需要用消元法: 第一个方程加上第二个方程乘以 $2$, 得到:\n",
    "$$\n",
    "5 x = - 4 \\Rightarrow x = -0.8.\n",
    "$$\n",
    "然后代入第二个方程, 得到\n",
    "$$\n",
    "y = 1.1.\n",
    "$$\n",
    "这是大家初中就会的算法. 现在的问题是, 如果问题推广到 $n$ 元一次**线性**方程组: (1.1.1), 那么我们学过的消元法能否跟着推广从而得到求解 $n$ 元一次**线性**方程组的解? 这里强调线性, 是关于未知数, 我们只有加法(减法)和数乘运算, 不能出现 $x^2$, $xy$ 或者 $\\sin(x)$ 这样的非线性运算."
   ]
  },
  {
   "cell_type": "markdown",
   "id": "701b2ed6",
   "metadata": {},
   "source": [
    "要回答这个问题, 我们首先应该解决一个技术问题: 我们如何书写, 记录一个 $n$ 元一次线性方程组? 要简单, 又要直观, 还要便于演算. "
   ]
  },
  {
   "cell_type": "markdown",
   "id": "9ef5017d",
   "metadata": {},
   "source": [
    "显然, 这个问题无关智商, 我们需要的是一套符号系统和约定. 本书的约定是一致的, 见 (1.1.1) 后对各符号的规定. 但是大家要重视这些符号和约定. 对每一个符号, 我们必须搞清楚它的实际含义, 并且要能够和具体的问题建立联系. 要特别避免两种情况: \n",
    "\n",
    "1. 觉得某个记号不重要, 直接跳过； 当你跳过一个你不理解的数学记号的时候, 你就失去了正确理解记号之后全部内容的可能性；切记! 切忌!\n",
    "2. 仗着自己惊人的记忆力强行将定义背下. 这就意味这你数学学得越多, 你记忆的负担就越重! 最后你不得不顶着一个巨大的记忆数据库去做逻辑推理. 这是人工智能的玩法, 而不是人类的玩法. 事实上先进的人工智能都已经开始摒弃这种模型. 正确的做法是通过建立和正确的具体问题的联系, 正确理解和接受这个符号. 当你看到一个符号时, 能自然联系起和它有关的信息和问题. 这样才会数学学得越多, 人也就越强.\n",
    "\n",
    "\n",
    "\n",
    "大家总是觉得数学很难, 特别是当技术文献中出现一行数学公式的时候, 真就是一种我和 paper 只能活一个的感觉. 这里最关键的 这里提示一下, 如果一个一元二次方程组是齐次的, 那么几何上看, 它的解是什么情况? 这个直观能推广到 $n$ 元么? \n",
    "\n",
    "要完整回答第一个问题, 我们还缺少一些必要的基础知识, 留待以后再讨论. 我们这里先接受课本对 $n$ 元一次线性方程组的解(以后用线性方程组或方程组简称)的规定. 提示: 在二元情形, 举一个方程组不相容的例子? 这种直观能推广到 $n$ 元么? "
   ]
  },
  {
   "cell_type": "markdown",
   "id": "7806a4f9",
   "metadata": {},
   "source": [
    "当一本数学书开始像我们的课本第一页(太残暴了)这样做定义疯狂输出时. 我们必须慢慢地推敲每一个定义的实际意义. 这里起首\"设 $\\mathbb{P}$ 是一个**数域**\", 我们可以扫边上的二维码来了解什么是数域. 简单地说, 数域就是对加减乘除封闭的数集(忽略被零除这个例外). 我们常见的数域有**有理数**($\\mathbb{Q}$), **实数**($\\mathbb{R}$)和**复数**($\\mathbb{C}$). **整数**不是数域. 在本书中, 我们默认数域为 $\\mathbb{R}$, 但因为抽象为一般数域 $\\mathbb{P}$, 所以我们大部分的结论, 对 $\\mathbb{Q}$, $\\mathbb{R}$ 和 $\\mathbb{C}$ 都成立. 这就是抽象的意义. 这里大家也可以看到, 数学的定义和符号, 蕴含非常丰富的信息. 对于第一次接触到的内容, 务必细细品味和推敲, 切忌求快. 这里慢就是快! "
   ]
  },
  {
   "cell_type": "markdown",
   "id": "abfbfccf",
   "metadata": {},
   "source": [
    "比如: 什么是**齐次线性方程组**? 写一个看看? 照之前的套路, 先来一个最具体最简单的. 二元一次齐次线性方程组:\n",
    "$$\n",
    "\\left\\{\n",
    "\\begin{array}{rcl}\n",
    "3 x + 4 y &=& 0 \\\\\n",
    "x - 2 y &=& 0.\n",
    "\\end{array}\n",
    "\\right.\n",
    "$$\n",
    "它的几何直观显然是: 两条过原点的直线. 而这个方程的意义是: 它们有公共点么? 废话! 至少有原点!"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "0dd20bf2",
   "metadata": {},
   "outputs": [],
   "source": [
    "x = np.linspace(-2,2,100)\n",
    "fig, ax = plt.subplots()\n",
    "line1, = plt.plot(x, (- 3*x) / 4, label = r'$3 x + 4 y = 0$')\n",
    "line2, = plt.plot(x, (- x) / -2, label = r'$x - 2 y = 0$')\n",
    "plt.plot([0.0], [0.0], 'o')\n",
    "ax.legend(handles=[line1, line2])"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "55595a8d",
   "metadata": {},
   "source": [
    "还能有更多的\"**解**\"么? 这里课本上强调解是一个**解集**的概念. 也就是说, 在某些情况下, 解可能不止一个. 那么在\"二元一次(这是线性的一个特例)齐次方程组\"的情况下, 可能有更多的解么? 我们从几何意义不难想到, 两条线可以重合? 或者它们实际上是同一条线(但方程**系数**不同)? 比如:\n",
    "$$\n",
    "\\left\\{\n",
    "\\begin{array}{rcl}\n",
    "3 x + 4 y &=& 0 \\\\\n",
    "6 x + 8 y &=& 0.\n",
    "\\end{array}\n",
    "\\right.\n",
    "$$\n",
    "形式上不同, 但几何上其实就是同一条线. 所以这样的线性方程组, 实际上有无穷多解: 解集为\n",
    "$$\n",
    "\\left\\{x \\left| x = -\\frac{4}{3} y, y \\in \\mathbb{R} \\right.\\right\\}.\n",
    "$$\n",
    "实际上就是这条重合的线上的每一点. 而在表示解集的时候, 我们用\"自由的\" $y$ 来表示 $x$, 这种方式便于理解 $x$ 和 $y$ 之间的约束关系.  \n",
    "\n",
    "现在, 我们在二元一次方程组的情形下有了一个具体的结论: 齐次线性方程组, 要么只有零解($x = y = 0$); 要么有无穷多解. 这个结论能推广到一般 $n$ 元么? "
   ]
  },
  {
   "cell_type": "markdown",
   "id": "a4449917",
   "metadata": {},
   "source": [
    "这里我们遇到了面对 $n$ 元或者 $n$ 维问题的一个最严重的问题: 我们画不出来了. 或者说, 我们失去了二维几何意义那样简单, 直接的直观概念. 这一点也体现了代数的必要性. 我们需要一个手段, 在失去这种直接直观的意义下, 能继续正确展开我们的推理和思考. 这种手段就是形式代数. 但我们不是说就彻底放弃直观了, 而是通过代数的形式推理, 我们试图建立一种超越我们的认知和常识的全新的直观. 这才是代数学习的核心意义, 也是当我们掌握了代数学这一思维工具, 能够体现出我们个体优势的根本原因. "
   ]
  },
  {
   "cell_type": "markdown",
   "id": "b734be57",
   "metadata": {},
   "source": [
    "那怎么才能发展出 $n$ 元的理解呢? 我们先猥琐发育一把. 考虑一下我们的二元理论能推导到三元么? 因为三维图形, 我们还是能勉强想象一下的. 大家都应该学过解析几何了吧? 我们考虑一个三元一次的齐次线性方程组是什么:\n",
    "$$\n",
    "A x + B y + C z = 0.\n",
    "$$\n",
    "我故意写成这个样子, 而不是我们书上的\n",
    "$$\n",
    "a_{11} x_1 + a_{12} x_2 + a_{13} x_3 = 0.\n",
    "$$\n",
    "因为前者学过解析几何的同学就知道, 这是一个平面啊. 过原点的平面. 因为显然 $x = y = z = 0$ 是它的解. 登时, 我们对我们结论的一部分: 对 $n$ 元问题, 原点是它的一个解充满了信心. 我们这时突然发现, 其实根本不用画图, 任何 $n$ 元齐次线性方程:\n",
    "$$\n",
    "a_{11} x_1 + a_{12} x_2 + \\cdots + a_{1n} x_n = 0,\n",
    "$$\n",
    "显然都满足零解: $x_1 = x_2 = \\cdots = x_n = 0$. 这一点其实从观察方程形式就能看到. 但是, 什么时候它会有无穷多解? 从二元的情形, 我们看到, 无穷多解的情况, 是它看上去有两个方程, 实际几何上只有一个方程. 这里我们先指出, 之前的讨论其实有问题. 我们一直到现在, 总是默认: 二元一次线性方程组, 必然由两个方程构成. 也就是如果用课本的语言, 我们认为 $m = n = 2$. 然而课本并没有规定 $m$ 和 $n$ 必须相等. 比如当二元一次方程组时, $n = 2$, 但 $m$ 可以是 $1$. 也就是看起来和实际上, 都只有一个方程. 这种情况下, 当然可以有无穷多解. 那么 $m > n = 2$ 呢? 这种鬼畜的情形, 要看实际上是不是真的有两个以上的独立方程. 如果有, 那么当然也只有零解. 实际上也可能出现, 看上去 $m > n$, 但实际上只有一个独立方程(可以考虑一下). 这样还是可能有无穷多解(一定包含零解). 而在三维, 如果出现 $3$ 个以上真正独立的方程, 那么会出现多个平面交于原点的情况, 这就是只有零解的情况. 而如果实际上只有不到 $3$ 个独立方程, 比如只有两个. 那么两个平面会交出一条过原点的直线, 这直线上的全部点都是解, 因此是无穷多个. "
   ]
  },
  {
   "cell_type": "markdown",
   "id": "057fd9fc",
   "metadata": {},
   "source": [
    "现在我们可以整理一下我们对 $n$ 元结论的猜测了: $n$ 元齐次线性方程组, 要么只有零解(实际上的独立方程个数大于等于 $n$), 或者有包含零解在内的无穷多解(实际上的独立方程个数小于 $n$). 这个结论, 大家是不是觉得怪怪的? 因为\"实际上的独立方程个数\"这个不太像个正经的数学定义. 出现这种情况, 就说明我们的定义和符号, 不够用了. 我们需要仔细推敲, 究竟这个概念是什么意思? 把它明确无误地用现有的符号系统表达出来, 可能需要做适当的扩展. 然后我们还需要在新体系下证明这一结论. 所以, 到这里, 我们先停下我们的讨论, 因为这些讨论已经涉及到好几周以后的内容了. 我们只是借用这一概念的展开, 领略一下, 数学符号的丰富内涵, 如何将它们和具体的问题建立联系, 如何用抽象的逻辑推理严格我们的结论, 再将其扩展到更一般的具体问题上. 这整个过程, 就是代数学. 这是一个动态的, 灵活的, 生机勃勃的过程, 不要只是背背背. 最后强调一下, 齐次线性方程组确实是一个容易被疏忽的知识点, 但如果你完整地经历了刚才的讨论, 那么相信你已经在你自己的知识体系中建立了一个\"什么是齐次线性方程组\"的直观, 这个才是我们课程学习, 讨论的最终目的: 建立数学概念的正确直观! 然后你还需要通过反复的推敲和作业练习, 来确保你的直观是正确的! 这需要付出一定的努力练习, 过程未必很快乐, 但我保证结果一定会让你在面对代数没学好的同学时充满优越感."
   ]
  },
  {
   "cell_type": "markdown",
   "id": "3ec0bb37",
   "metadata": {},
   "source": [
    "以上讨论了书本的一些概念, 但真正重要的内容是我们应该如何学习数学! 顺便说一句, 既然有了解集的概念, 那么当解集为空时, 就是课本所定义的**不相容**."
   ]
  },
  {
   "cell_type": "markdown",
   "id": "cb07ff3f",
   "metadata": {},
   "source": [
    "现在我们来继续讨论, 如何将消元和代入的方法, 推广到 $n$ 元一次的一般情形."
   ]
  },
  {
   "cell_type": "markdown",
   "id": "ce59d783",
   "metadata": {},
   "source": [
    "## 1.2 同解变形与阶梯线性方程组\n",
    "\n",
    "我们注意到, 消元法的本质在于, 将一个形式上一般的方程组, 化简为一个形式简单的方程组, 也就是在一些方程组, 自变量**依次**减少, 最后出现只有一个自变量的方程, 然后从它开始求解, 得到结果之后再代入其它方程, 依次得到全部子变量的值(或范围关系, 参考前面只有一个实质方程的形式). 在化简的过程中, 关键是不能改变方程组的**解集**. 也就是要做**同解变形**."
   ]
  },
  {
   "cell_type": "markdown",
   "id": "0429d504",
   "metadata": {},
   "source": [
    "感谢历史上 Gauss 等猛人, 他们提供了课本上描述的三种基本变形手段: 互换, 倍乘和倍加. 这三种手段加在一起, 称为线性方程组的**初等变换**. 这里我们稍微理解一下, 互换两个方程(或者互换两个未知量的项, 这个在下注中)本质上没做任何事情, 只是书写需要. 因此解不变是显然的. 倍乘, 从我们前面的经验, 两个系数是倍乘关系的方程实质上是同一个方程. 或者说如果:\n",
    "$$\n",
    "x_1 = c_1, x_2 = c_2, \\cdots, x_n = c_n\n",
    "$$\n",
    "是方程\n",
    "$$\n",
    "a_{i1} x_1 + a_{i2} x_2 + \\cdots + a_{in} x_n = b_i\n",
    "$$\n",
    "的解, 那么它必然也是\n",
    "$$\n",
    "c \\cdot a_{i1} x_1 + c \\cdot a_{i2} x_2 + \\cdots + c \\cdot a_{in} x_n = c \\cdot b_i\n",
    "$$\n",
    "的解. 没毛病啊. 而倍加只要考虑如果\n",
    "$$\n",
    "x_1 = c_1, x_2 = c_2, \\cdots, x_n = c_n\n",
    "$$\n",
    "也是\n",
    "$$\n",
    "a_{j1} x_1 + a_{j2} x_2 + \\cdots + a_{jn} x_n = b_j\n",
    "$$\n",
    "的解, 那么它必然是方程\n",
    "$$\n",
    "(a_{i1} + a_{j1}) x_1 + (a_{i2} + a_{j2}) x_2 + \\cdots + (a_{in} + a_{jn}) x_n = (b_i + b_j)\n",
    "$$\n",
    "的解. 于是我们就验证了引理 1. 现在我们就有了三种方程同解变形的手段, 可以拿来做消元. "
   ]
  },
  {
   "cell_type": "markdown",
   "id": "8213742a",
   "metadata": {},
   "source": [
    "看课本例 1-(1), 我们的目的是通过消元化简方程. 首先看看哪个方程放第一行. 这里课本建议我们放一个 $a_{11} = 1$ 的, 因为这样手算简单. 这里你可以选择第一行或第二行都行. 接下去我们考虑如何有次序地消元, 不要一阵乱消. 比如, 我们总是以 $x_1$, $x_2$, $\\cdots$ 的次序依次消元. (提示: 这个次序可以改变, 但那就意味着互换两个未知量的项.)以这个次序, 我们先用第一个方程, 消全部其余方程中的 $x_1$ 项. 对第二行, 直接减去第一个方程, 这是初等变换的特例, 课本用 $R_2 - R_1$ 记录, 表示 \"row 2 - row 1\". 对第三行, 减去第一个方程乘以 $3$, 记作 $R_3 - 3 R_1$. 然后 $R_4 - 5 R_1$. 这样完成对 $x_1$ 的消元. "
   ]
  },
  {
   "cell_type": "markdown",
   "id": "0837f12f",
   "metadata": {},
   "source": [
    "这个例子很特殊, 消去 $x_1$ 的时候, 顺便把这三个方程中的 $x_2$ 也消去了. 这种好事不常有. 这样我们接下去就只要消 $x_3$ 了. 这里注意, 第一个方程, 我们不要再去动它. 不然来来回回没完没了. 实际上, 我们目前是打算构造一个算法, 用于一般性的消元. 第一个方程的信息我们已经充分利用了, 所以现在就把它放在那里. 我们接下去考虑继续用第二个(新)方程, 它一定不含 $x_1$(在这个很特殊的例子中, 连 $x_2$也没有了), 然后我们可以放心地用它去消剩余方程中的 $x_3$. 继续用之前的记号, $R4 - 2 R_3$, 和 $R_5 - R_2$.  "
   ]
  },
  {
   "cell_type": "markdown",
   "id": "2ce18190",
   "metadata": {},
   "source": [
    "课本的最后一步初等变换 $-\\frac{1}{7} R_2$ 和 $R_1 - 5 R_2$ 只是一种形式整理, 并不涉及消元. 最后得到一种鬼畜的情形, 最后两行 $0 = 0$ 很对, 但似乎也完全没有用的, 因为它们没有提供任何有意义的信息. 这让我们联系之前的讨论, 到底什么是 \"独立的方程个数\" 啊? 这里我们看到, 在同解的初等变换下, 有两行方程完全失去约束的意义. 因此我们可以认为, 这里只有**两个独立的有约束意义的方程**. 而从剩下的第一行和第二行我们可以看到, 这个方程组是应该是一个解集. "
   ]
  },
  {
   "cell_type": "markdown",
   "id": "1bf15363",
   "metadata": {},
   "source": [
    "再看例 1-(2). 它的不同之处在于最后得到的形式中, 居然出现 $0 = 1$. 这个显然是不可能的. 因此告诉我们这个方程组其实**不相容**, 无界. 但是这一行在方程约束上, 和第四行 $0 = 0$ 一样无意义. 因此这里真正有意义的方程组个数, 还是两个."
   ]
  },
  {
   "cell_type": "markdown",
   "id": "ad84482e",
   "metadata": {},
   "source": [
    "现在我们来考虑, 如何用初等变换, 对一个一般的 $n$ 元一次方程组做消元. 得到更一般的结论. 我们接受课本 (1.1.1) 的符号体系, 线性方程组写作:\n",
    "$$\n",
    "\\left\\{\n",
    "\\begin{array}{rcl}\n",
    "a_{11} x_1 + a_{12} x_2 + \\cdots + a_{1n} x_n &=& b_1, \\\\\n",
    "a_{21} x_1 + a_{22} x_2 + \\cdots + a_{2n} x_n &=& b_2, \\\\\n",
    "\\cdots &\\cdots& \\cdots \\\\\n",
    "a_{m1} x_1 + a_{m2} x_2 + \\cdots + a_{mn} x_n &=& b_m. \\\\\n",
    "\\end{array}\n",
    "\\right.\n",
    "$$"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "269f9aab",
   "metadata": {},
   "source": [
    "现在考虑一般意义下的消元. 我们可以利用初等变换, 先用第一个方程, 消去其余各方程中的 $x_1$ 项. 具体来说, 对第 $i$ 个方程:\n",
    "$$\n",
    "a_{i1} x_1 + a_{i2} x_2 + \\cdots + a_{in} x_n = b_i, \n",
    "$$\n",
    "我们用第一个方程乘上 $-\\frac{a_{i1}}{a_{11}}$ 再加到第 $i$ 个方程上(倍加), (即 $R_i - \\frac{a_{i_1}}{a_{11}} R_1$, $i = 2, 3, \\cdots, m$)得到:\n",
    "$$\n",
    "0 \\cdot x_1 + \\left(a_{i2} - a_{12} \\cdot \\frac{a_{i1}}{a_{11}}\\right) x_{2} \n",
    "+ \\left(a_{i3} - a_{13} \\cdot \\frac{a_{i1}}{a_{11}}\\right) x_{3} + \\cdots \n",
    "+ \\left(a_{in} - a_{1n} \\cdot \\frac{a_{i1}}{a_{11}}\\right) x_{n} \n",
    "= b_i - b_1 \\cdot \\frac{a_{i1}}{a_{11}}\n",
    "$$\n",
    "其中, $i = 2, 3, \\cdots, m$. 这里注意我们要假设 $a_{11} \\neq 0$. (如果等于零咋办? 换一个啊! 课本从收算方便的角度, 建议大家换一个 $a_{11}$ 尽量简单, 比如 $a_{11} = 1$, 但如果以后你要写程序的话, 实际上这个$a_{11}$ 应该尽可能的大). 注意这个过程是完全机械的. 对这个规则, 可以用数学的语言, 或者用程序来描述. "
   ]
  },
  {
   "cell_type": "markdown",
   "id": "7810db48",
   "metadata": {},
   "source": [
    "其实为了简化书写, 我们可以将全部的方程系数集中起来, 写成:\n",
    "$$\n",
    "A = \\left(\n",
    "\\begin{array}{cccc}\n",
    "a_{11}&a_{12}&\\cdots&a_{1n}\\\\\n",
    "a_{21}&a_{22}&\\cdots&a_{2n}\\\\\n",
    "\\cdots&\\cdots&\\cdots&\\cdots\\\\\n",
    "a_{m1}&a_{m2}&\\cdots&a_{mn}\n",
    "\\end{array}\n",
    "\\right).\n",
    "$$\n",
    "在今后, 这个形式会被称为**矩阵**. 而同样的, 为方便起见, 将全部右端项也集中成一个向量的形式:\n",
    "$$\n",
    "b = (b_1, b_2, \\cdots, b_m)^T.\n",
    "$$\n",
    "一般来说, $b$ 应该是一个列向量, 但我们用一个 $T$ 表示转置运算, 主要是为了省书写空间. 比如课本例1, 用矩阵和向量的记号, 可以写成:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "226f297a",
   "metadata": {},
   "outputs": [],
   "source": [
    "A = np.array([[1, 2, 5, -1], [1, 2, -2, 1], [3, 6, 8, -1], [5, 10, 11, -1]], dtype='double')\n",
    "b = np.array([-1, 3, 1, 3], dtype='double')"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "2bd1c72a",
   "metadata": {},
   "source": [
    "程序和数学具有一样的逻辑. 它们可以描述同一个规则. 这里, 我们可以将上述矩阵和向量, 自动转化为线性方程组的形式:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "ac2a3764",
   "metadata": {},
   "outputs": [],
   "source": [
    "def print_Equations(A, b):\n",
    "    m = len(A)\n",
    "    n = len(A[0])\n",
    "    for i in range(m):\n",
    "        for j in range(n):\n",
    "            if (np.abs(A[i][j]) < eps):\n",
    "                print(\"            \\t\", end='')\n",
    "            else:\n",
    "                print(\"%+f x%d\\t\" % (A[i][j], j + 1), end='')\n",
    "        print(' =\\t%+f\\n' % b[i], '\\n')        "
   ]
  },
  {
   "cell_type": "markdown",
   "id": "4fe9c9c8",
   "metadata": {},
   "source": [
    "用这个程序规则, 我们不难在程序中输出我们的方程组:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "c5c3b0c2",
   "metadata": {
    "scrolled": true
   },
   "outputs": [],
   "source": [
    "print_Equations(A, b)"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "316bbd19",
   "metadata": {},
   "source": [
    "这样我们甚至可以把消元第一步写成程序:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "deb933bc",
   "metadata": {},
   "outputs": [],
   "source": [
    "def Gaussian_Eliminate(A, b):\n",
    "    m = len(A)\n",
    "    n = len(A[0])\n",
    "    for i in range(1, m):\n",
    "        for j in range(1, n):\n",
    "            A[i][j] = A[i][j] - A[0][j] * A[i][0] / A[0][0]\n",
    "        b[i] = b[i] - b[0] * A[i][0] / A[0][0]\n",
    "        A[i][0] = 0\n",
    "    return (A, b)"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "238cc477",
   "metadata": {},
   "source": [
    "我们看这个规则在具体问题例 1 下的效果:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "d588e5b9",
   "metadata": {},
   "outputs": [],
   "source": [
    "(A, b) = Gaussian_Eliminate(A, b)\n",
    "print_Equations(A, b)"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "3a430c69",
   "metadata": {},
   "source": [
    "注意我们的规则是一般的, 因此完全可以对一个随机产生的线性方程组套用此规则:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "f120ea89",
   "metadata": {},
   "outputs": [],
   "source": [
    "A = [[np.random.random() for col in range(4)] for row in range(4)]\n",
    "b = np.random.random(4)\n",
    "print_Equations(A, b)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "3e038f58",
   "metadata": {
    "scrolled": true
   },
   "outputs": [],
   "source": [
    "A1 = np.array(A)\n",
    "b1 = np.array(b)\n",
    "(A1, b1) = Gaussian_Eliminate(A1, b1)\n",
    "print_Equations(A1, b1)"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "3cbce362",
   "metadata": {},
   "source": [
    "在消 $x_1$ 时, 有可能把 $x_2$ 甚至更多的未知量消掉, 因此消完 $x_1$ 可能接着是消 $x_{j_2}$, 这个 $j_2$ 可能是 $2$, $3$, $\\cdots$, 所以书上把实际消元的次序, 用 $x_{j_1}$, $x_{j_2}$, $\\cdots$, $x_{j_r}$ 表示. 并规定\n",
    "$$\n",
    "1 = j_1 < j_2 < \\cdots < j_r \\leq n.\n",
    "$$"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "62fca138",
   "metadata": {},
   "source": [
    "现在假设我们已经对 $x_{j_i}$ 完成了消元, 那么现在方程应该长这个样:\n",
    "$$\n",
    "\\begin{array}{rrllclclcl}\n",
    "a^{(0)}_{11} x_1 + \\cdots + &a^{(0)}_{1j_2} x_{j_2} + \\cdots + &a^{(0)}_{1j_i} x_{j_i} + \\cdots + &a^{(0)}_{1j_{i + 1}} x_{j_{i + 1}} &+ \\cdots + & a^{(0)}_{1j_{i + 2}} x_{j_{i + 2}} &+ \\cdots + &a^{(0)}_{1n} x_n & = & b_1^{(0)}, \\\\\n",
    "&a^{(1)}_{2j_2} x_{j_2} + \\cdots + &a^{(1)}_{2j_i} x_{j_i}  + \\cdots + &a^{(1)}_{2j_{i + 1}} x_{j_{i + 1}} &+ \\cdots + & a^{(1)}_{2j_{i + 2}} x_{j_{i + 2}} &+ \\cdots + &a^{(1)}_{2n} x_n & = & b_2^{(1)}, \\\\\n",
    "&\\ddots &\\vdots&\\vdots&\\cdots&\\vdots&\\cdots&\\vdots&=&\\vdots\\\\\n",
    "&&a^{(i - 1)}_{ij_i} x_{j_i}  +\\cdots + &a^{(i - 1)}_{ij_{i + 1}} x_{j_{i + 1}} &+ \\cdots + & a^{(i - 1)}_{ij_{i + 2}} x_{j_{i + 2}} &+ \\cdots + &a^{(i - 1)}_{in} x_n & = & b_i^{(i - 1)} \\\\\n",
    "&&&a^{(i)}_{i+1j_{i+1}} x_{j_{i + 1}} &+ \\cdots + & a^{(i)}_{i + 1j_{i + 2}} x_{j_{i + 2}} &+ \\cdots + &a^{(i)}_{i+1n} x_n & = & b_{i + 1}^{(i)} \\\\\n",
    "&&&a^{(i)}_{i+2j_{i+1}} x_{j_{i + 1}} &+ \\cdots + & a^{(i)}_{i + 2j_{i + 2}} x_{j_{i + 2}} &+ \\cdots + &a^{(i)}_{i+2n} x_n & = & b_{i + 2}^{(i)} \\\\\n",
    "&&&\\vdots&\\cdots& \\vdots&\\cdots& \\vdots & = & \\vdots \\\\\n",
    "&&&a^{(i)}_{mj_{i + 1}} x_{j_{i + 1}} &+ \\cdots + & a^{(i)}_{mj_{i + 2}} x_{j_{i + 2}} & + \\cdots + &a^{(i)}_{mn} x_n & = & b_m^{(i)} \n",
    "\\end{array}\n",
    "$$\n",
    "这里上标 $^{(k)}$ 表示经过 $k$ 步消元的改变. $k = 0$ 表示没有改变过. 要达到上述形式, 经过了 $i$ 步消元."
   ]
  },
  {
   "cell_type": "markdown",
   "id": "617c4283",
   "metadata": {},
   "source": [
    "现在我们可以利用初等变换, 用第 $i + 1$ 行方程, 消去第 $i + 2 \\sim m$ 行各方程中的 $x_{j_{i + 1}}$ 项. 具体来说, 我们用第 $i + 1$ 个方程乘上 $-\\frac{a^{(i)}_{k~j_{i + 1}}}{a^{(i)}_{i+1~j_{i + 1}}}$ 再加到第 $k$ 个方程上(倍加, $R_k -\\frac{a^{(i)}_{k~j_{i + 1}}}{a^{(i)}_{i+1~j_{i + 1}}} R_{i + 1}$), $k = i + 2, i + 3, \\cdots, m$, 得到:\n",
    "$$\n",
    "\\begin{array}{rrllclclcl}\n",
    "a^{(0)}_{11} x_1 + \\cdots + &a^{(0)}_{1j_2} x_{j_2} + \\cdots + &a^{(0)}_{1j_i} x_{j_i} + \\cdots + &a^{(0)}_{1j_{i + 1}} x_{j_{i + 1}} &+ \\cdots + & a^{(0)}_{1j_{i + 2}} x_{j_{i + 2}} &+ \\cdots + &a^{(0)}_{1n} x_n & = & b_1^{(0)}, \\\\\n",
    "&a^{(1)}_{2j_2} x_{j_2} + \\cdots + &a^{(1)}_{2j_i} x_{j_i}  + \\cdots + &a^{(1)}_{2j_{i + 1}} x_{j_{i + 1}} &+ \\cdots + & a^{(1)}_{2j_{i + 2}} x_{j_{i + 2}} &+ \\cdots + &a^{(1)}_{2n} x_n & = & b_2^{(1)}, \\\\\n",
    "&\\ddots &\\vdots&\\vdots&\\cdots&\\vdots&\\cdots&\\vdots&=&\\vdots\\\\\n",
    "&&a^{(i - 1)}_{ij_i} x_{j_i}  +\\cdots + &a^{(i - 1)}_{ij_{i + 1}} x_{j_{i + 1}} &+ \\cdots + & a^{(i - 1)}_{ij_{i + 2}} x_{j_{i + 2}} &+ \\cdots + &a^{(i - 1)}_{in} x_n & = & b_i^{(i - 1)} \\\\\n",
    "&&&a^{(i)}_{i+1j_{i+1}} x_{j_{i + 1}} &+ \\cdots + & a^{(i)}_{i + 1j_{i + 2}} x_{j_{i + 2}} &+ \\cdots + &a^{(i)}_{i+1n} x_n & = & b_{i + 1}^{(i)} \\\\\n",
    "&&&&& a^{(i + 1)}_{i + 2j_{i + 2}} x_{j_{i + 2}} &+ \\cdots + &a^{(i + 1)}_{i+1n} x_n & = & b_{i + 2}^{(i + 1)} \\\\\n",
    "&&&&& \\vdots&\\cdots& \\vdots & = & \\vdots \\\\\n",
    "&&&&& a^{(i + 1)}_{mj_{i + 2}} x_{j_{i + 2}} & + \\cdots + &a^{(i + 1)}_{mn} x_n & = & b_m^{(i + 1)} \n",
    "\\end{array}\n",
    "$$"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "e8781564",
   "metadata": {},
   "source": [
    "这是一个算法过程, 这个方法被称为 Gauss 消元法. 对于一般的方程, 我们如果假设, 每次消元只能消掉一个未知量, 那么可以将这个过程写成如下的程序:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "2a5ef93a",
   "metadata": {},
   "outputs": [],
   "source": [
    "def Gaussian_Eliminate(A, b):\n",
    "    m = len(A)\n",
    "    n = len(A[0])\n",
    "    for k in range(m - 1):\n",
    "        for i in range(k + 1, m):\n",
    "            for j in range(k + 1, n):\n",
    "                A[i][j] = A[i][j] - A[k][j] * A[i][k] / A[k][k]\n",
    "            b[i] = b[i] - b[k] * A[i][k] / A[k][k]\n",
    "            A[i][k] = 0\n",
    "    return (A, b)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "0096916e",
   "metadata": {},
   "outputs": [],
   "source": [
    "(A1, b1) = Gaussian_Eliminate(A1, b1)\n",
    "print_Equations(A1, b1)"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "e117b7d4",
   "metadata": {},
   "source": [
    "对于书上更一般的结论, 我们不能期待总是会出现\"上三角型\"消元矩阵, 只能保证出现一个**阶梯型线性方程组**. 且这个方程组和原始线性方程组是同解的. 这个就是课本的定理1结论. 我们以上的分析过程, 实际上就是对定理 1 的一个严格证明. 课本结论中的 (1.2.1), 有两个信息值得关注, 首先是**阶梯头**, 也就是实质性方程一共有 $r$ 个. 现在我们终于对什么是\"独立的实质性方程个数\"有了一个构造性定义: 如果通过初等变换进行 Gauss 消去, 得到的阶梯头个数 $r$ 就是实质性独立方程的个数. 或者说, 在阶梯型线性方程组中去掉那些 $0 = 0$ 或 $0 = 1$ 等鬼畜方程以后剩下的方程数量. 第二个值得注意的信息是第 $r + 1$ 行方程, 如果有, 必然是 $0 = 1$ 或 $0 = 0$ 型. 前者, 意味着方程组无解, 而后者, 方程组有无穷多解. 那么是否可能出现方程组有且仅有一个解? 直观上似乎这个才是应用中最常见的情况? "
   ]
  },
  {
   "cell_type": "markdown",
   "id": "e3bcbae5",
   "metadata": {},
   "source": [
    "这种情况就是 $r = n$ 的情形! 比如我们上面的随机方程组例子. 它的第 $r$ 行就是最后一行, 而且直接给出了 $x_r$ 的值. 然后我们可以代回倒数第二个方程, 解得 $x_3$, 如此重复, 直到得到全部解. 更一般地, 如果对 $n$ 个方程构成的 $n$ 元一次方程组, 且在消去过程中, 所有的 $a_{kk} \\neq 0$, $k = 1, 2, \\cdots, n$, 则我们可以得到:\n",
    "$$\n",
    "x_i = \\left(b_i - \\sum_{j = i + 1}^n a_{ij} x_j\\right) / a_{ii}, \n",
    "$$\n",
    "这里 $i = n, n - 1, \\cdots, 1$. 这个过程因此形象地被称为\"回代过程\"."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "d827f9cb",
   "metadata": {},
   "outputs": [],
   "source": [
    "def back_substitute(A, b):\n",
    "    m = len(A)\n",
    "    n = len(A[0])\n",
    "    x = np.zeros(n)\n",
    "    for i in range(m - 1, -1, -1):\n",
    "        s = 0\n",
    "        for j in range(i, n):\n",
    "            s += x[j] * A[i][j]\n",
    "        x[i] = (b[i] - s) / A[i][i]\n",
    "    return x    "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "c6bce834",
   "metadata": {
    "scrolled": true
   },
   "outputs": [],
   "source": [
    "x = back_substitute(A1, b1)\n",
    "print(x)"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "8021cda8",
   "metadata": {},
   "source": [
    "这里有必要写一个程序验证一下解是否正确(或者手算验证): "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "f7104714",
   "metadata": {},
   "outputs": [],
   "source": [
    "def check_solution(A, b, x):\n",
    "    m = len(A)\n",
    "    n = len(A[0])\n",
    "    error = b\n",
    "    res = 0\n",
    "    for i in range(m):\n",
    "        for j in range(n):\n",
    "            error[i] -= A[i][j] * x[j]\n",
    "        res += np.abs(error[i])\n",
    "    return res"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "dda25da1",
   "metadata": {},
   "outputs": [],
   "source": [
    "check_solution(A1, b1, x)"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "6bf31273",
   "metadata": {},
   "source": [
    "这里再验证一下, 所得 $x$ 是否为原方程组的解:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "80697eef",
   "metadata": {},
   "outputs": [],
   "source": [
    "check_solution(A, b, x)"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "1ab8d156",
   "metadata": {},
   "source": [
    "那么对于更一般的情形, 也就是 $r < n$ 的情形, 比如例 1-(1)? 注意我们的程序是解不了这个问题的."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "f1a76048",
   "metadata": {},
   "outputs": [],
   "source": [
    "A = np.array([[1, 2, 5, -1], [1, 2, -2, 1], [3, 6, 8, -1], [5, 10, 11, -1]], dtype='double')\n",
    "b = np.array([-1, 3, 1, 3], dtype='double')\n",
    "Gaussian_Eliminate(A, b)"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "b06174cb",
   "metadata": {},
   "source": [
    "课本的引理 2回答了这个问题. 大家可以看到, 我们上面的回代法, 其实是引理 2 结论的一个特例. 如果 $r = n$ 成立, 总是能得到唯一解. 但是当 $r < n$ 等条件下时, 实际上我们会得到无穷多的解. 因为此时, $x_{j_{r + 1}}, x_{j_{r + 2}}, \\cdots, x_{j_n}$ 这 $n - r - 1$ 个未知量, 事实上是**没有约束**的. 因此, 我们叫它们**自由未知量**. 它们事实上可以取数域 $\\mathbb{P}$ 中的任意数. 一旦自由未知量取定, 则其它未知量, 通过 (1.3.7) 是唯一确定的. 但自由未知量的任意取法, 导致原方程组实际上存在无穷多解. 对以上结论的系统性整理, 就是我们的课本定理 2."
   ]
  },
  {
   "cell_type": "markdown",
   "id": "17a480af",
   "metadata": {},
   "source": [
    "例 2 给出了例 1-(1)的最终回代求解. 注意! 真正的独立变量, 是阶梯头! 比如这里, 我们不能把 $x_2$ 定为独立变量. 它必须是自由变量! 因为, 在消元时 $x_2$ 和 $x_1$ 一起被消去, 蕴含了 $x_2$ 和 $x_1$ 其实是相关的, 这个概念我们会在后续的课程中展开讨论. 这里我们要明确引理 2 中, 独立变量是 $x_{j_1}$, $x_{j_2}$, $\\cdots$, $x_{j_r}$ 这些引理 1 中的**阶梯头**变量, 而剩下的变量 $x_{j_{r + 1}}$, $x_{j_{r + 2}}$, $\\cdots$, $x_{j_n}$, 是**自由变量**, 可以在数域 $\\mathbb{P}$ 中随意取. "
   ]
  },
  {
   "cell_type": "markdown",
   "id": "bf3a0054",
   "metadata": {},
   "source": [
    "## 1.4 矩阵及其初等变换\n",
    "\n",
    "我们之前已经引出了矩阵的记号. 用矩阵来表达线性方程组的消去过程, 非常方便. 我们首先观察, 在矩阵形式下, 初等变换的定义仍然相容. 注意到我们的消去法, 实际上是要影响到右端向量 $b = (b_1, b_2, \\cdots, b_m)$ 的, 因此, 将系数矩阵 $A$ 和右端向量 $b$ 一起写成一个分块的增广矩阵, 非常方便:\n",
    "$$\n",
    "(A | b) = \\left(\n",
    "\\begin{array}{cccccc}\n",
    "a_{11}&a_{12}&\\cdots&a_{1n}&|&b_1\\\\\n",
    "a_{21}&a_{22}&\\cdots&a_{2n}&|&b_2\\\\\n",
    "\\cdots&\\cdots&\\cdots&\\cdots&|&\\cdots\\\\\n",
    "a_{m1}&a_{m2}&\\cdots&a_{mn}&|&b_m\n",
    "\\end{array}\n",
    "\\right).\n",
    "$$\n",
    "我们可以在这个增广矩阵上记录和演算高斯消元的过程. 参考课本例 3 和例 4. 这两题的风格和考试接近."
   ]
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 3",
   "language": "python",
   "name": "python3"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 3
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython3",
   "version": "3.7.9"
  },
  "widgets": {
   "application/vnd.jupyter.widget-state+json": {
    "state": {},
    "version_major": 2,
    "version_minor": 0
   }
  }
 },
 "nbformat": 4,
 "nbformat_minor": 5
}
